#coding:utf-8

import random
import numpy as np
import pylab as p

def teste(a,b,c,x,y,z,r1,r2):
	return pow((x-a),2)+pow((y-b),2)+pow((z-c),2)-pow((r1+r2),2)

def volume(v):
	max = [-100,-100,-100]
	min = [100,100,100]	
	
	for i in range(0,n):
		ind=3*i
		if v[ind]+ri[i] > max[0]:
			max[0] = v[ind]+ri[i]
		if v[ind+1]+ri[i] > max[1]:
			max[1] = v[ind+1]+ri[i]
		if v[ind+2]+ri[i] > max[2]:
			max[2] = v[ind+2]+ri[i]			
			
		if v[ind]-ri[i] < min[0]:
			min[0] = v[ind]-ri[i]
		if v[ind+1]-ri[i] < min[1]:
			min[1] = v[ind+1]-ri[i]
		if v[ind+2]-ri[i] < min[2]:
			min[2] = v[ind+2]-ri[i]
		
	return ( ( ( abs(max[0]-min[0]) )*( abs(max[1]-min[1]) )*( abs(max[2]-min[2]) ) ),min,max )
	
def gen_dir(nvar):
	dirs = np.zeros( (nvar+1,nvar) )
	for i in range(0,nvar):
		dirs[i,i] = 1
		
	dirs[nvar,:] = [-1]*nvar
	
	return dirs

def gen_esferas(n,max_r,max_pt):
	xi = [ max_pt*random.random() for j in range(0,n) ]
	yi = [ max_pt*random.random() for j in range(0,n) ]
	zi = [ max_pt*random.random() for j in range(0,n) ]
	
	ri = [ int( max_r*random.random() ) for j in range(0,n) ]
	
	return (xi,yi,zi,ri)

n=10

res = gen_esferas(n,15,50)
xi = res[0]
yi = res[1]
zi = res[2]
ri = res[3]

class Util:
	MAX_IT = 5000
	FATOR_1 = 0
	FATOR_2 = 1.5

def soma_v(v1,v2):
	res = [0]*len(v1)
	for i in range(0,len(v1)):
		res[i] = v1[i]+v2[i]
	return res
	
def interno_v(v1,v2):
	res = 0
	for i in range(0,len(v1)):
		res += v1[i]*v2[i]
		
	return res	
	
def escalar_v(v1,k):
	res = [0]*len(v1)
	for i in range(0,len(v1)):
		res[i] = v1[i]*k
		
	return res	
	
def norma(v):
	res = 0
	for i in range(0,len(v)):
		res+=pow(v[i],2)
	res = np.sqrt(res)
	return res

	
def Ei(i,v):
	ind_i=3*i
	res = 0
	for k in range(0,n):
		ind_k=3*k			
		
		if i==k:
			continue
			
		parc = pow( (v[ind_i]-v[ind_k]),2 ) + pow( (v[ind_i+1]-v[ind_k+1]),2 ) + pow( (v[ind_i+2]-v[ind_k+2]),2 ) - pow( (ri[i]+ri[k]),2 )
		res+=pow(min(parc,0),2)
		
	return res
	
def Fi(i,v):
	ind =3*i
	res = Util.FATOR_1*(pow(v[ind],2)+pow(v[ind+1],2)+pow(v[ind+2],2))
	res += Util.FATOR_2*Ei(i,v)
	
	return res
	
def F(v):
	res = 0
	
	max = [-100,-100,-100]
	min = [100,100,100]	
	for i in range(0,n):
		res+=Fi(i,v)
		
		ind = 3*i
		if v[ind]+ri[i] > max[0]:
			max[0] = v[ind]+ri[i]
		if v[ind+1]+ri[i] > max[1]:
			max[1] = v[ind+1]+ri[i]
		if v[ind+2]+ri[i] > max[2]:
			max[2] = v[ind+2]+ri[i]			
			
		if v[ind]-ri[i] < min[0]:
			min[0] = v[ind]-ri[i]
		if v[ind+1]-ri[i] < min[1]:
			min[1] = v[ind+1]-ri[i]
		if v[ind+2]-ri[i] < min[2]:
			min[2] = v[ind+2]-ri[i]
		
	res += ( ( abs(max[0]-min[0]) )*( abs(max[1]-min[1]) )*( abs(max[2]-min[2]) ) )
		
	return res
	
def main():
	delta = 1.0
	epson = 1e-6
	it = 0
	
	nvar = 3*n
	dirs = gen_dir(nvar)
	num_dir = np.arange( 0,len(dirs) )
	x0 = [0]*nvar
	
	f = F(x0)
	
	while (delta>epson and it<=Util.MAX_IT):
		#print it,f
		found_dir = False
		mask = np.array( [True]*len(dirs) )
		
		while True:
			dirs_left = num_dir[mask]
			if len(dirs_left)==0:
				break
				
			ndir_choice = random.choice( dirs_left )
			mask[ndir_choice] = False
			
			d = dirs[ndir_choice]
			
			xk = soma_v(x0,escalar_v(d,delta))
			f_ = F( xk )		
			
			if f_ < f:
				x0 = xk
				f=f_
				found_dir=True			
				break		
				
		if not(found_dir):
			delta/=1.2
			
		it+=1
	
	print it,delta
	return x0

def print_scilab(x):
	print "v=[",
	for i in range(0,len(x)):
		if i%3==0 and i!=0:
			print ";",
		print x[i],
	print "]"

	print "ri=",ri

	resp = volume(x)
	print "caixa=[ ",resp[1],"],[",resp[2],"]"
	print "Volume=",resp[0]	
	
def iniciar():
	x = main()
	#x=[ 51.5646972656, -19.0, 5.0 , -49.2451171875, -67.0,78.625 , -0.644897460938, 0.625, 0.0 ]
	print_scilab(x)

def varias_vezes():
	min = 1000000
	xr = None
	for i in range(0,50):
		x = main()
		v = volume(x)[0]
		if v < min:
			min = v
			xr = x
			
	print_scilab(xr)
	
varias_vezes()